Mott phases and superfluid-insulator transition of dipolar spin-three bosons in an 

optical lattice: implications for 52 Cr atoms 
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We study the Mott phases and superfluid-insulator transition of spin-three bosons in an optical 
lattice with an anisotropic two dimensional optical trap. We chart out the phase diagrams for 
Mott states with n = 1 and n = 2 atoms per lattice site. It is shown that the long-range dipolar 
interaction stabilizes a state where the chains of the ferromagnetically aligned spins run along the 
longer trap direction while the spin ordering is staggered between nearby chains, leading to an 
antiferromagnetic ordering along the shorter trap direction. We also obtain the mean-field phase 
boundary for the superfluid-insulator transition in these systems and study the nature of spin 
ordering in the superfluid state near the transition. We show that, inside the superfluid phase and 
near the superfluid-insulator phase boundary, the system undergoes a first order antiferromagnetic- 
ferromagnetic spin ordering transition. We discuss implications of our results for 52 Cr atoms and 
suggest possible experiments to detect different phases in such systems. 

PACS numbers: 



I. INTRODUCTION 

Several experiments on ultracold trapped atomic gases 
have opened a new window onto a plethora of intrigu- 
ing phases of quantum matter A gas of bosonic atoms 
in an optical or magnetic trap has been reversibly tuned 
between superfluid (SF) and insulating ground states by 
varying the strength of a periodic potential produced by 
standing optical waves. This transition has been ex- 
plained on the basis of the Bose-Hubbard model with 
on-site repulsive interactions and hopping between near- 
est neighboring sites of the lattice^. As long as the atom- 
atom interactions are small compared to the hopping am- 
plitude, the ground state remains superfluid. In the op- 
posite limit of a strong lattice potential, the interaction 
energy dominates and the ground state is a Mott insu- 
lator (MI) when the density is commensurate, with an 
integer number of atoms localized at each lattice site. 

In the above-mentioned experiments, the spins of the 
atoms were neglected either due to freezing of the spin 
degrees of freedom in the presence of a magnetic trap 
or due to the use of spinless atoms. Later, use of op- 
tical trap for confining atoms led to realization of spin- 
one and spin-two bosonic systems where the spins of the 
atoms are dynamical degrees of freedom^ 5 -. Furthermore, 
possibilities of creating two-component bosonic mixture, 
which is equivalent to bosons with pseudospin half, has 
also been discussed 6 -. In all of these systems, the bosonic 
atoms are shown to have spin-dependent contact inter- 
action which can lead to variety of spin ordered Mott 
ground states 6 -^^ ' 10 : 11 . The superfluid-insulator transi- 
tions in these systems has also been studied^ ' 11 : 12 . 



More recently, a Bose-Einstein condensate of 52 Cr 
atoms has been realized^ 3 -. These Cr atoms, in contrast 
to their alkali counterparts, have the electronic configura- 



tion [Ar] 3d 5 4s 1 which leads, for the isotopes with no nu- 
clear spin, to a net spin of 3 (in units of K) and to a mag- 
netic moment of 6^b, where /is is the Bohr magneton 1 " 4 *. 
Consequently, these atoms experience a much stronger 
long-range dipolar interaction compared to their alkali 
counterparts. Besides, again in stark contrast to other 
bosonic alkali atoms, in a Cr condensate, the scatter- 
ing lengths in different angular momentum channel are 
quite differem—. Since the Cr atoms have spin 3, the 
relevant scattering lengths (as) occurs in angular mo- 
mentum channels S = 0, 2, 4, 6. All the as except ao are 
determined and have the values a 2 = — 7as, 04 = 58a b 
and a§ = 112as, where as is the Bohr radiu a 14 ' 15 . Such 
different scattering lengths lead to much stronger spin- 
dependent contact interaction potentials in Cr compared 
to its alkali counterparts 15 . In fact, in the presence of 
an optical trap, where the spins become dynamical de- 
grees of freedoms, the spin-dependent contact interac- 
tions overwhelm the long-ranged dipolar interaction en- 
ergy, and lead to various spin-nematic phase s 15 : 16 . In 
particular, it has been predicted in Ref. [l5| that biaxial 
nematic phases may be realized in Cr condensates^. 

In this work, we study the Mott phases and superfluid- 
insulator transition of spin-three bosons in an optical lat- 
tice with an anisotropic two dimensional optical trap. 
The main results reported in this paper are as follows. 
First, we map out the phase diagrams for Mott states 
with n = 1 and n — 2 atoms per lattice site. In par- 
ticular, we demonstrate that both the contact and the 
dipolar interactions are important for obtaining the cor- 
rect Mott ground states. It is shown that the presence of 
a long-range dipolar interaction and an anisotropic two- 
dimensional (2D) trap lead to a Mott ground state where 
the chains of the ferromagnetically (FM) aligned spins 
run along the longer trap direction while the spin order- 
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ing is staggered between nearby chains, leading to an an- 
tiferromagnetic (AFM) ordering along the shorter trap 
direction. Second, we obtain the superfluid-insulator 
transition (SIT) phase boundary of this system using 
mean-field theory and discuss the spin ordering of the 
superfluid near and away from the transition. We demon- 
strate that, in contrast to the alkali bosonic atoms with 
weak spin-dependent contact interaction and negligible 
dipolar interactions, the spin-three bosons with dipolar 
interaction exhibit a first order AFM-FM transition for 
the spin ordering in the condensate as one moves away 
from the SIT boundary into the superfluid phases. Fi- 
nally, we suggest future experiments on 52 Cr to test our 
theoretical predictions. 

We start with the following Hamiltonian 



H = H' K + Ti-'o + H' D 
H' k '■ 



dr V>t(r 
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where i/>t (r) is the boson creation operator with spin pro- 
jection a = {—3,-2,-1,0,1,2,3} at position r, V(r) is 
the optical lattice potential, S™ b denotes the element of 
the three spin-3 rotation matrices, sum over all repeated 
indices are implied, and u(r — r') is the two-body inter- 
atomic potential given by 



u(r 



S(r - r')(g P + g 2 P 2 + g 4 P 4 + g 6 P 6 ), (5) 



where g s = ^h 2 a s /M and P s = EL-s \S,m)(S, m\ is 
the projection operator for the Bosons. Using the iden- 
tities 



1 = PV + P2 + P4 + P6, 
(S-S') = -12P - 9P 2 - 2P 4 + 9P 6 , 
(S-S') 2 = 144P + 81P 2 +4P 4 + 81P 6 , 

we can rewrite the contact interaction as 



(6) 



u(r-r') = S(r-r') [u + + u 2 (S • S') + u 3 (S • S') 2 ] , 

(7) 

where u = (~llg 2 + 8I.94 + 7g 6 )/77, u x = g Q + (-55g 2 + 
27 5 4 -5 ff6 )/33, u 2 = ( 56 - 52 )/18, u 3 = {g 2 /126-g 4 /77 + 
96/198). The coefficient of the dipolar term 114 is known 
to be small (i.e. u 4 <C uq, u\,u 2 , u 3 ) in 52 Cr atomic. 

For ultracold atoms in an optical lattice, it is well- 
known that the energy eigenstates are Bloch wave func- 
tions and a superposition of these Bloch states yields a 
set of Wannier functions that are well localized on the 
individual lattice sites for deep enough lattices^ 8 . Fur- 
ther, in these systems both the kinetic and interaction 



energies involved in the dynamics are small compared to 
the excitation energies to higher single particle bands in 
the limit of deep optical lattices necessary to realize the 
Mott insulating statei*i&. We shall therefore expand the 
boson field operators in the Wannier basis keeping only 
the lowest single particle band, V'a(r) = Yli b i a w(r — ri), 
where bi a is the spin-3 boson annihilation operator at 
site i with spin projection a. After a few straightfor- 
ward manipulations, we then obtain the effective single 
band extended Bose-Hubbard Hamiltonian on a lattice 
Hb = Hr + Ho + Hd, where 



Hk = 



Ho = 



(ij) a 

^5>(rW) + ^£(S 2 -12 fk) 

i i 



(8) 
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Hd = 



i k,l—{x,y,z} 
^4 / Sj • Sj 

o qv qv ( V i ~ ( r i ~ r j)v 



(9) 



(10) 



Here hi = 

Si = Eat 



Ea b la b ia 



bLSabbi 



J2ab b ia( S kSl)abhb 



is the boson density at site i, 
is the spin operator and Ni^i = 
is the 'nematic' operator. The 
hopping coefficient t is given by t — J dr w*(r — 

ri) ^2l7^ 2 + ^( r )) w ( r — r i) an( i the are gi ven 
by Ui = Uif dr|w(r)| 4 , for i = {0,1,2,3}, and U 4 = 
u 4 J dr\w(r)\ 2 Jdr'\w(r')\ 2 . 

The Hamiltonian Hb is the starting point of our study 
in the subsequent sections. In what follows, we shall 
obtain the ground states of Hb — M J2i 

function of 

Ui/Uq, U 2 /Uq, U3/U0, fJ,/Uo, and t/Uo, where // denotes 
the chemical potential for a fixed weak dipolar interaction 
U 4 /U <C 1, f^j={i,2,3}/^o- Note that since the value of 
ao and hence U\ is not experimentally determined, the 
position of 52 Cr in this generalized phase diagram is not 
known. However, we expect that some of the general 
features obtained from this phase diagram will apply to 
52 Cr. Also, we would like to point out that in principle 
we expect the chemical potential fj, to be space-dependent 
due to the presence of the trapping potential. In this 
work, we have ignored this effect of the optical trap which 
is a standard approximation in the literature for large 
trapsi£. 

The rest of the paper is organized as follows. In Sec. [Hi 
we study various possible Mott phases of the Hb- This 
is followed, in Sec. IIII1 by a mean-field study of the SIT 
in this model. Then in Sec. IIV1 we suggest experiments 
to test our predictions. Finally, we conclude in Sec. |Vj 
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II. MOTT PHASES 

In this section, we obtain various possible Mott phases 
when there are integer n particles per site. In the Mott 
limit, the kinetic energy of atoms are negligible and we 
shall take t = in the rest of this section. We shall 
concern ourselves only with the case of n — 1 and n = 2 
atoms per site in this work. The calculations can be 
generalized, in principle, to higher n values. The phase 
diagram for the on-site terms is sketched in Sec. Ill Al while 
the effect of the dipolar interaction is presented in Sec. 

EU 

A. On-site Interaction 

In order to obtain the Mott phases, we use the fact 
that the dipolar interaction term Hd (Eq. IT0|) is small 
compared to the on-site interaction terms in Ho (Eq. 
IH]). Thus, we first neglect the dipole term Hd and find 
out the ground-state of Ho by minimizing the variational 
energy 

E v {n) = (^ v (n)\n \^v(n)), (11) 

where the on-site variational wavefunction \ty v (n)} for n 
bosons per lattice site is given by 

s 

= [] ^ c »,(™,s,m)|^;S',m} l ,(12) 
i S=S to tai m=— S 

where Stotai denotes all the possible total spin that can be 
obtained by adding spins of n spin-three bosons and m is 
the azimuthal quantum number corresponding to S. For 
example for n = 2, the possible values of Stotai = 0,2,4,6. 

Details of the calculation of E v (n) are sketched in Ap- 
pendix [A] Here we present the phase diagram obtained 
by numerical minimization of E v (n) for n = 1 and n = 2 
as a function of U3/U0 and U2/U0 for several representa- 
tive values of fi/Uo and for Ui/Uq = ±0.2. These values 
of /i and U\ are chosen for clarity; we have checked that 
there are no other phases that occur for other values of 
H and U\. In order to find the lower boundary in U3/U0 
for the n — 2 states, we also minimize E v (n = 3). How- 
ever, since we are mostly concerned with the n = 1 and 
71 = 2 states and also for clarity, we have not indicated 
the different S regions of n = 3 in the phase diagrams 
shown in Figs. [U f2J [3J and |U From our calculations, we 
find that Ho do not lift the azimuthal degeneracy of the 
atoms and the ground states obtained from this on-site 
variational wavefunction are degenerate. The only excep- 
tion to this occurs for the singlet ground state (S = 0) 
that exists in a narrow region of the phase diagram. For 
all other phases, the on-site interactions determine only 
the S value of the ground state, thus leaving a 25 + 1 
degeneracy corresponding to different possible m values. 
This is where the dipole interaction comes into play. The 
presence of the dipolar interaction fixes the direction of 



l-i/Uo = 0.01, Ui/Uo = 0.2 
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FIG. 1: Mott phase diagram obtained by minimization of 
E v (n) for n = 1 and n = 2 as a function of U2/U0 and U3/U0 
for n = 0.01(7o and U1/U0 = 0.2. Note that the azimuthal 
degeneracy is not lifted and hence the ground states for non- 
zero S have a 2S + 1 fold degeneracy. 
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FIG. 2: Same as in Fig. □ for y, = O.8U0 and Ui/Uq = 0.2 



the spins and gives us a unique ground state. Notice 
that the weakness of the dipolar interactions compared 
to the on-site term Ho ensures that the former would 
not change the value of S in the ground state. 



H/Uo = I.OMJUq = 0.2 
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n = 2 S = 4 


n = 2 S = 6 


. „ = I =~ 


n = 3 
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FIG. 3: Same as in Fig. [T] for fx = U and U1/U0 = 0.2 
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li/Vo = l-O.UJUo = -0.2 



n = 1 S = 3 










/ n = 2 S = 4 


n = 2 S = 6 y 




n = 2^ S =_2___ 


~~> n = 2 5' = 




n = 3 





-0.05 0.0 0.05 0.1 0.15 

FIG. 4: Same as in Fig. U for [X = U and [/i/C/ = -0.2 



B. Dipolar Interaction 

To investigate the effect of the dipolar interaction, we 
now restrict ourselves to the lowest S manifold which has 
been determined by minimizing E v (n). However, even 
within this manifold, we need to construct a variational 
wavefunction to obtain the ground state of TLd- Since the 
dipole interaction is long-ranged, such a trial wavefunc- 
tion need not be a simple on-site product wavefunction. 
To guess the form of this wavefunction, we first look at 
the dipolar interaction for classical spins in an optical lat- 
tice and find out the ground state for such spins. Next, 
we choose a variational wavefunction which is compatible 
with the symmetry of this classical ground state. 

In the classical limit, the spins can be parameterized 
as Si = S (sin 9i cos 4>i,siii6i sin fa, cos 0i). Also, due to 
the form of the dipolar interaction (Eq. flQ|) and the two- 
dimensional nature of our system, it is easy to see that 
the ground state, in the classical limit, shall have all spins 
aligned in the XY plane. Thus in the ground state, each 
classical spin can be represented by a single angle cf>. This 
simplifies the computational procedure enormously and 
enables us to study the ground states of classical spins for 
different trap configurations with both open and periodic 
boundary conditions. 

To find the minimum energy configuration for the clas- 
sical spins, we minimize 'Hoi'pi) f° r both square and 
rectangular trap geometries. For square traps with open 
boundary conditions or rectangular traps with periodic 
boundary conditions in both x and y directions, the bulk 
region exhibits a swirling pattern described by four inter- 
penetrating sublattices (see FigJ5]) with the spins making 
an angle /3 with the lattice at each site (see Fig.|B]). More- 
over, simulations with periodic boundary conditions in 
both trap directions show that the system energy is de- 
generate for all values of (3. Consequently, for square 
traps, dipolar interactions are not sufficient to lift the 
spin degeneracy. However, this situation is difficult to 
realize in realistic finite size systems with optical traps 
where a rectangular trap with open boundary condition 
is a more realistic choice. For the rest of this work, we 
shall therefore concentrate on such geometries. In this 
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FIG. 5: 14 x 14 lattice with the open boundary condition. No- 
tice that the spins in the bulk region are following a "swirling" 
pattern that can be described by four interpenetrating sub- 
lattices. 
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FIG. 6: Typical spin pattern for classical spins with dipolar 
interactions in the square lattice with the periodic boundary 
condition. The ground state has infinite degeneracy corre- 
sponding to different values of /3. 



case, the spins of the bulk region tend to align with the 
longest trap direction as shown in Fig(7] Consequently, 
in the classical sense, the spin degeneracy is completely 
lifted by dipolar interactions in rectangular geometries. 
The ground state for the classical spins for rectangular 
trap with open boundary condition thus have an anti- 
ferromagnetic(AFM) order along the shorter trap direc- 
tion. Note that such an order corresponds to an inter- 
penetrating four-sublattice structure shown in Fig.UJwith 
[3 = ±7r/2 when x is taken as the longest trap direction. 

The lesson we learn from the study of the dipolar 
Hamiltonian in the classical limit is that for rectangular 
traps with open boundary conditions, the ground state 
has four interpenetrating sublattices as seen in Fig. [7] 
Therefore, a reasonable minimal trial wavefunction for 
finding the ground state of Tip in the quantum limit must 
at least have a four-sublattice structure. With this obser- 
vation, we choose the following variational wavefunction 
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FIG. 7: 14 x 7 lattice with the open boundary condition. 
Notice that the spins in the bulk region are aligned along the 
longer trap direction. This pattern can be described via four 
interpenetrating sublattices. 



for n = 1 

|*dipolo(« = 1)) 



n n x ( i3 ) 

A={A,B,C,D] ieA 
3 

Cj,(n=l,S=3,m)|l;3,m)i, 

m=— 3 



where {A, B, C, D} are the four sublattices. Similarly for 
n = 2, one gets 



|*di P ole(n = 2)) = 



n n x ( m ) 

A={A,B,C,D} ieA 
S 

X] X! c i,{n=2,S,m)\%S,m)i. 
S={0,2,4,6}m=-S 



We 



then 



minimize E^ lpolc (n) = 
(*dipole(n)| Hd |*di P ole(n)) to obtain the ground 
state for the system. We have carried out this numerical 
minimization for rectangular lattice sizes up to 16 x 8 
with periodic boundary conditions in one of the two 
directions. For all values of 5 ^ and for all lattice 
sizes up to 16 x 8, we find again that the spins are 
aligned parallel to the longest trap direction and have 
the same four interpenetrating sublattice structure as 
shown in Fig. [7] Since this ground state turns out to be 
generic for all 5 7^ 0, we expect to find such a ground 
state for 52 Cr atoms irrespective of its position in the 
generalized phase diagram for n = 1. For n = 2, the 
ground state will either be the four interpenetrating 
sublattice structure mentioned previously (if 5 ^ 0) or 
a singlet (5=0). 



III. SUPERFLUID-INSULATOR TRANSITION 
A. Phase Boundary 

In this section, we obtain the superfluid-insulator 
phase boundary starting from the n = 1 Mott phase 



within a mean-field theory. To allow density fluctuations 
in the superfluid state, we construct a variational wave- 
function that is a superposition of n = 0, n = 1 and n = 2 
Mott states. Note that the variational wavefunction here 
needs to have the four sublattice structure mentioned in 
Sec. lIIBl to correctly capture the n = 1 Mott state. Such 
a wavefunction is given by 



i*>o = n ii x 

A={A,B,C,D} ieA 

(ci,(„=o)|rc = 0), +\n= l)i + \n = 2)i), 



(15) 



where 



\n = l) l = C *,(«=l,S=3,m)&!,JO} (16) 



771 = — 3 



and 



\n = 2)i = X X! c i,(n=2,s,m)\2;S,m)i 

S'={0,2,4,6} m=-S 
S 

— X X c i,{n=2,S,m) 

S={0,2,4,6} m=-S 



(/ b 



where ti{s, a +b=m) are Clebsch-Gordan coefficients. The 
normalization of |^)o leads to the condition 



c i,(™=0)| 2 + J]] c i,(n=l,S=3,m) I 
m— — 3 
S 

+ X X l c i,(n=2,S : m)| 2 - 

S={0,2,4,6} m=-S 



(18) 



Using this wavefunction, the expectation value of the 
kinetic term is given by 

E K = o(*|«x|*)o = -t ]T 2 5R(A* Q A ja ), (19) 

(i,j),a 

where the superfluid order parameter A* a = o^l&ij^o 
can be expressed in terms of coefficients of the variational 
wavefunction >F)o as 

A* a = c i,(n=l,S=3,a) C i,(n=0) 

S 3 

+ 2J X! X! C 7*(n=2,S,rn') 

S={0,2,4,6} m'=—S be m=-3 

xc i,(n=l,S=3,m)^(S,b+c=m')(SbaScm + S bm S C a)- 

(20) 

Since the interaction terms TLo and Jiu (Eqs. [5] and 
ITUl) do not mix states with different n, the total varia- 
tional energy can be written as 



E 



Mott 



'n = 1) 



SE h + SE n 



(21) 
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where SE p and SEh which are the energy costs of adding 
a particle and hole respectively to the Mott state, are 
given by 



SE n = - 



(*di P oic(n = 2)| 
(Ho +H d - /i^^j|* d ip i c (n = 2)) 

i 

(^dipolcO = 1)1 

7-to + T-Ld - M ^ j l*di P oic("- = 1)) 

i 

(*dipolc(™ = 1)1 



(22) 



N 



Ho + T~Ld - /Lt^njj |^di P oic(^ = 1)) 



SE h 



(23) 

The mean-field ground state of the system can now be ob- 
tained, as a function oft/Uo and /j/Uq for representative 
values of Ui/U , U2/U0, U 3 /U and U^/Uq, by numerical 
minimization of E var . Note that as t exceeds a critical 
value t c , the energy gain from Ek (Eq. [19]) exceeds the 
energy cost of adding a particle (SE p ) or a hole (5E h ) to 
the Mott state. At this point, the ground state of the 
system occurs with non-zero superfluid order parameter 

We obtain t c by minimizing E var (Eg. EH]) . In doing so, 
to reduce computational time, we only keep the two low- 
est manifolds of S states in the Mott phase corresponding 
to |n = 2)i (Eq. [T"7]) . The justification of this procedure 
comes from the following. At the superfluid-insulator 
phase boundary, the hopping coefficient t reaches a crit- 
ical value for which it becomes energetically favorable to 
add (or remove) a particle to the n = 1 Mott state. When 
such a particle is added, the resultant state has n = 2 
and can, in principle, have superposition of S — 0, 2, 4, 6 
states due to the presence of the dipole term. However, 
when the dipole interaction is weak compared to the on- 
site interaction, such superposition becomes energetically 
costly. Thus the resultant lowest-energy n = 2 state is 
almost entirely localized within the two lowest S mani- 
folds for t ~ t c . We have numerically checked for a few 
cases keeping all the S states that this is indeed the case. 

The mean-field superfluid-insulator phase boundary, so 
obtained, is shown in Fig. [8] for representative param- 
eter values Ui = Q, U 2 = 0MU o , U 3 = 0.02[/ and 
U4 = 0.0125C/o- Here the two lowest n = 2 states cor- 
responds to S = 2 and S — 4 manifolds and we have 
kept these states for our minimization calculation. We 
note that compared to the usual Bose-Hubbard model 3 , 
the transition occurs for a relatively higher t c /U. This 
phenomenon can be understood to be the consequence 
of the spin-ordering of the Mott state which originates 
from the dipole interaction. For sake of clarity, we shall 
consider the situation where the Mott state is destabi- 
lized by addition of holes. The other case, where the 
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FIG. 8: The mean-field phase boundary for the superfluid- 
insulator transition with the n — 1 Mott state in an optical 
lattice. Here U 2 = 0.05C/ , U 3 = 0.02(7 and U 4 = 0.0125l/" . 
Notice the relatively higher values of t c /Uo compared to the 
usual cases, which is a consequence of spin order in the Mott 
phase (see text). 



Mott state is destabilized by addition of particles, has a 
similar explanation. 

Consider the Mott state of the system with the spin or- 
der as shown in Fig. [7] It turns out that the coefficients 
Ci corresponding to spins in nearby chains have opposite 
signs because the spin directions are staggered from one 
chain to the other. Since A* a ~ c?(n=i,S=3,a) c i,(«=o)> 
the sign of varies between different sites. Conse- 
quently, the sign of A*Aj in Eq. [19] can vary so that 
the kinetic energy is 'frustrated'. This leads to a reduc- 
tion of the kinetic energy and hence one needs a higher 
t value to destabilize the Mott phase. We note that this 
phenomenon is generic and is a consequence of dipolar 
interaction between the spins, which stabilizes the above- 
mentioned four-sublattice order. 



B. AFM-FM transition 

In this section, we study the spin-ordering of the su- 
perfluid phase near the transition. The strategy for such 
a study is again minimization of E val followed by com- 
putation of (S x ) and (S y ) using the variational ground 
state wave function. However, even before such a com- 
putation is carried out, one can guess the result from a 
simple physical picture, which we now present. 

We have already noted in Sec. lIII Al that the AFM spin- 
order along the shorter trap direction, which minimizes 
the dipole energy in the Mott state, is in competition 
with maximization of kinetic energy near the superfluid- 
insulator transition. Thus as we go deeper in the super- 
fluid phase, we expect that there would be a transition 
from a state with such a spin-order to a FM state where 
all the spins point towards x (the longer trap direction). 
Such a transition from the translational-symmetry bro- 
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ken AFM state to a uniform FM state is expected to be 
first order. The exact location of this transition depends 
on the relative strength of t, ^ and U4. As shown in 
Fig.EJ the AFM-FM transition follows closely the MI-SF 
transition for values of h/Uq smaller than (/i/^o)t max an d 
departs from the MI-SF transition boundary for larger 
values of h/Uq. This departure from the MFSF bound- 
ary can be explained by noticing that the Mott state is 
destabilized primarily by addition of holes (particles) for 
fi/Uo close to 0(1). For intermediate values of //, the 
destabilization of the Mott phase occurs, within mean- 
field theory, with finite amplitudes for addition of both 
holes and particles. Consequently, for a fixed n/Uo close 
to 0, as we increase t/Uo, the average number of bosons 
per site decreases and (5 X ) decreases accordingly. There- 
fore, the effective dipole interaction becomes weaker and 
taSm, the value of t for which the kinetic energy gain wins 
over the optimization of the dipolar energy, stays close to 
t c . However, for a fixed n/Uo close to 1, the average num- 
ber of bosons per site increases as t is increased beyond 
t c and (S x ) also slowly increases. Thus, for [i > /i ( m«, 
the dipole interaction is robust to an increase in t/Uo so 
that i a f m becomes larger and deviates significantly from 
the Mott-superfluid phase boundary. 

We have also computed the variation of £ a fm with C/4, 
as shown in Fig. [10] for representative values of parame- 
ters U 3 = 0.02[/o, U 2 = 4U4 and n = 0.05U a . We find 
that i a f m approaches t c for small values of U±/Uq. The 
AFM-FM transition occurs simultaneously with the MI- 
SF transition for U4 < 0.00125C/o, leading to a first order 
MI-SF transition in this regime. Such dependence of t a f m 
on U& can be understood by noting that a reduction of C/4 
directly decreases the dipolar energy. As a result, when 
U4 is reduced, the AFM spin ordering becomes less ro- 
bust to increases in t because it becomes easier to take 
advantage of the kinetic energy gain at the expense of 
the dipolar energy. Thus £ a f m decreases, and ultimately 
merges to t c for U 4 /U Q < 0.00125. 

Finally, we investigate the fate of the spin-order as 
we go deeper in the SF region. We note at the outset 
that since the approximation involved in our minimiza- 
tion scheme loses its accuracy as we go deeper into the 
superfluid phase, we do not expect our analysis to be 
quantitatively accurate in this regime, but only expect 
some qualitative features to be captured. To this end, we 
plot (Sx) as a function of t/Uo with \i/Uo = 0.05 (with all 
other parameters being the same as in Fig. [8]) as shown in 
Fig. [H] Note that here t c ~ 0.16C/ and i afm ~ 0.17J7 , 
so that t c and i a f m are quite close to each other. As we 
further increase t, we find that for t/Uo > 0.27, we have 
a transition to the paramagnetic phase where (S x ) = 0. 
We have also computed the expectation value of the ne- 
matic order parameters in the FM state. The results of 
these computation are summarized in Fig. 1121 We find 
that for t < 0.23C/o and t > 0.27J7o, we have a uniaxial 
nematic phase. For 0.23 < t/Uo < 0.27, in the region 
over which (S x ) decrease from its maximum value to 0, 
we find a biaxial nematic phase, similar to one found in 
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FIG. 9: The AFM-FM and the MI-SF transition boundaries 
for the n = 1 Mott state in an optical lattice as a function of 
H/U . Here U 2 = 0.05[/ , Us = 0.02[/ and U* = 0.0125£/o. 
Notice that the AFM-FM transition follows closely the MI- 
SF transition for values of fi/Uo smaller than (/x/[/o)t™ ax and 
departs from the MI-SF transition boundary for larger values 
of n/U . 
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FIG. 10: The AFM-FM and the MI-SF transition boundaries 
for the n = 1 Mott state in an optical lattice as a function 
of U4/U0 for Us = 0.02C/0, U2 = 4C/ 4 and fi = 0.05[/ . No- 
tice that for Ui = 0.00125(7 , the AFM-FM transition occurs 
simultaneously with the MI-SF transition. 



Ref. [l5j. However, we note that at the point where the 
biaxial nematicity sets in, we are already deep inside the 
SF phase and our approximations of retaining the two 
lowest manifold of states may be inaccurate. 



IV. EXPERIMENTS 

The traditional way of examining the existence of su- 
perfluidity in trapped boson systems is to switch off the 
trap, let the cloud of atoms expand freely and image 
the expanding cloud. The momentum distribution of 
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FIG. 11: Total magnetization in the x direction (M^ otal = 
(%2 i Sf)/N, where the sum is over all lattice sites) for the 
Mott and superfluid phases as a function of t/Uo- Here /i = 
0.05C/ , U 2 = 0.05C/ , U 3 = 0.02t/ , and U 4 = 0.0125t/ o . Note 
that the AFM-FM transition occurs for i afm ~ 0.17(7 > t c ~ 
0.16(7o and, as we increase t, the FM order gives away to a 
paramagnetic superfluid for t > 0.27t/o- 



all the different phases. To achieve this distinction, we 
need to use the well-known technique of passing the ex- 
panding cloud through a Stern Gerlach magnet with the 
magnetic field along y (shorter trap direction)^. Such a 
Stern-Gerlach magnet will split the expanding cloud into 
two clouds of equal proportion if the atoms in the trap 
were in the AFM phase. If the atoms were in the FM 
phase, no such splitting will occur. Thus, a momentum- 
distribution measurement, along with the Stern Gerlach 
measurement described above can distinguish between 
the Mott- AFM (delocalized momentum distribution and 
two clouds), superfluid- AFM (localized momentum dis- 
tribution and two clouds) and superfluid-FM phases (lo- 
calized momentum distribution and single cloud). How- 
ever, we note that such measurements do not uniquely de- 
termine the details of the spin configuration in the AFM 
state. The signature of the details of the spin configu- 
ration can be obtained experimentally by measuring the 
spatial noise correlations of the expanding clouds 2 ^. Fi- 
nally, the uniaxial and biaxial nematic orders can also be 
detected using an imaging technique suggested in Ref. 
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FIG. 12: AFM, FM, superfluid and nematic phases as a func- 
tion of t/Uo for fi = 0.05C/ , U 2 = 0.05t/ , U s = 0.02[/ and 
U 4 = 0.0125f/ . Here we have plotted (T,i=i,L x Sf)/{W X ), 
where the sum is over the N x number of sites in a given chain 
along the x direction (the longer trap direction). It turns out 
that its magnitude in the superfluid phase varies with /i and 
the U's. 



V. CONCLUSION 

In conclusion, we have studied the Mott insulating 
phases and superfluid-insulator phase transitions of S = 
3 bosons in an optical lattice. In particular, we have con- 
sidered dipolar interaction between the bosons and the 
effect of rectangular optical trap. We have investigated 
possible Mott phases of such systems for integer n = 1 
and n = 2 atoms per lattice site and shown that the 
presence of dipolar interaction leads to a state with fer- 
romagnetic chains along the longer trap direction which 
are aligned antiparallel to each other, thus leading to an 
antifcrromagnetic order along the shorter trap direction. 
We have also presented a mean-field phase diagram for 
the superfluid-insulator transition in these systems and 
have shown that there is a generic first order AFM-FM 
transition in the superfluid region close to the superfluid- 
insulator phase boundary. We have also suggested real- 
istic experiments on 52 Cr to verify our predictions. 



the atoms inside the trap can then be inferred by look- 
ing at their position, or equivalently density, distribution 
in the expanded cloud. Since the momentum distribu- 
tion function of the atoms is characterized by the pres- 
ence/absence of coherence peaks in the superfluid/Mott 
insulating states, such a measurement serves as a qual- 
itative probe of the state of the atoms inside the trap 1 . 
In our proposed setup, however, such a simple expan- 
sion alone, which is not sensitive to the different spin 
states of atoms, will not be able to distinguish between 
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APPENDIX A: COMPUTATION OF THE while for n = 2 the energy is given by 

VARIATIONAL ENERGY FOR THE MOTT 

PHASES s 



E v {2) = ii U - 2 H) E \ C Un=2,S,m)\ 2 



E v (n) = <*„(n)|{ -^^^(Tii-l) 



S m=-S 

S 



We present in this section the expression of the varia- jj 2 
tional energy E v (n) (Eq. HIjl in greater detail. Following +— ^ Y \ c i,(n=2,S,m)\ 2 (S(S + 1) - 24) 

Eq. [9j the on-site energy is given by s m=—S 

+ U\ |Ci ! („ = 2,S=0,m=0)| 2 

tt ^ — ^(288 |c< ( n — 2 5 = o m =o) | 2 



162 ^ |c l! („ = 2,S=2,rn)| 2 



m=— 2 
4 



i ab q 

-Vl T AT 2 H -132rW] + 162 E K(«=2, S =6,m)| 2 )}- (A5) 

t k,l— {x,y,z} 

fi hi } . The Mott phases are then obtained by minimiz- 

(Al) m 6 E v {n) with respect to the variational parameters 

^i, (n,S,m) ■ 



To compute the on-site energy in the case of n = 1, we 
use the variational wavefunction 



3 

*»(!)> =11 E *,(n=i,s=3,tn}&tj0>, (A2) 

i rn — — 3 



while for n = 2, we use 



s 

1^(2)) = LI E E C *,(»=2,S,m) 

« S={0,2,4,6} m=-S 

(A3) 



where £(s,a+b=m) are Clebsh-Gordan coefheients. Using 
these wave functions, we find the energies for n = 1 and 
n = 2. For n = 1, the on-site energy is given by 



3 

^(1) = -A«E E k(n=l,S=3,m)| 2 > (A4) 
i m— — 3 



10 



M. Greiner, O. Mandel, T. Esslinger, T.W. Hansch and I. 
Bloch, Nature 415, 39 (2002). 

C. Orzel, A.K. Tuchman, M.L. Fenselau, M. Yasuda and 
M.A. Kasevich, Science 291, 2386 (2001). 
M.P.A. Fisher, P.B. Weichman, G. Grinstein, and D.S. 
Fisher, Phys. Rev. B 40, 546 (1989). 

J. Stenger, S. Inouye, D.M. Stamper-Kurn, H.-J. Miesner, 
A.P. Chikkatur, W. Ketterle, Nature 396, 345 (1998). 
A. Widera, F. Gerbier, S. Foiling, T. Gericke, O. Man- 
del and I. Bloch, Phys. Rev. Lett. 95, 190405 (2005); A. 
Widera, F. Gerbier, S. Foiling, T. Gericke, O. Mandel and 
I. Bloch, New J. Phys 8, 152 (2006). 

L.-M. Duan, E. Demler, and M. D. Lukin, Phys. Rev. Lett. 
91, 090402 (2003); E. Altman, W. Hofstetter, E. Demler, 
and M.D. Lukin, New J. Phys. 5, 113 (2003). 
A. Isacsson, M-C Cha, K. Senguptaand S.M. Girvin, Phys. 
Rev. B 72, 184507 (2005). 

F. Zhou and G. W. Semenoff, Phys. Rev. Lett. 97, 180411 
(2006). 

S. Powell and S. Sachdev, [cond-mat/0 608611 
R. Barnett, A. Turner and E. Demler, Phys. Rev. Lett. 
97, 180412 (2006); R. Barnett , A. Turner and E. Demler, 
|cond-mat/0611230| S.-K. Yip, |cond-mat/061117l"| 



1 A. Imambekov, M. Lukin, and E. Demler, Phys. Rev. A 
68, 063602 (2003); J.-S. Bernier, K. Sengupta and Y.B. 
Kim, Phys. Rev. B 74, 155124 (2006). 

12 J.-M. Hou and M.-L. Ge, Phys. Rev. A 67, 063607 (2003). 

13 A. Griesmaier, J. Werner, S. Hensler, J. Stuhler and T. 
Pfau, Phys. Rev. Lett 94, 160401 (2006). 

14 J. Stuhler, A. Griesmaier, T. Koch, M. Fattori, T. Pfau, 
S. Giovanazzi, P. Pedri and L. Santos, Phys. Rev. Lett 95, 
150406 (2005). 

15 R.B. Diener and T-L. Ho, Phys. Rev. Lett. 96, 190405 
(2006). 

16 L. Santos and T. Pfau, Phys. Rev. Lett. 96, 190404 (2006). 

17 Refs. [TH and differ regarding the presence of a biaxial 
nematic phase, but agree on the presence of more common 
uniaxial nematic phases in these systems. 

18 D. Jaksch, C. Bruder, J.I. Cirac, C.W. Gardiner and P. 
Zoller, Phys. Rev. Lett. 81, 3108 (1998). 

19 C. Schroll, F. Marquadt, and C. Bruder, Phys. Rev. A 70, 
053609 (2004). 

20 E. Altman, E. Demler, and M.D. Lukin, Phys. Rev. A 70, 
013603 (2004). 

21 I. Carusotto and E. Mueller, J. Phys. B 37, S115 (2004). 



